---
title: "Sustainable food systems"
author: "Harold Achicanoy"
output:
flexdashboard::flex_dashboard:
orientation: columns
social: menu
source_code: embed
---
```{r setup, include=FALSE}
suppressMessages(library(highcharter))
suppressMessages(library(dplyr))
suppressMessages(library(viridisLite))
suppressMessages(library(forecast))
suppressMessages(library(treemap))
suppressMessages(library(flexdashboard))
suppressMessages(library(VIM))
suppressMessages(library(plspm))
suppressMessages(library(FactoMineR))
suppressMessages(library(networkD3))
suppressMessages(library(igraph))
suppressMessages(library(plyr))
suppressMessages(library(cluster))
suppressMessages(library(plotly))
thm <-
hc_theme(
colors = c("#1a6ecc", "#434348", "#90ed7d"),
chart = list(
backgroundColor = "transparent",
style = list(fontFamily = "Source Sans Pro")
),
xAxis = list(
gridLineWidth = 1
)
)
```
Column {data-width=600}
-----------------------------------------------------------------------
### Sales Forecast
```{r}
AirPassengers %>%
forecast(level = 90) %>%
hchart() %>%
hc_add_theme(thm)
```
### Sustainable food system index by country
```{r}
data("worldgeojson")
indices <- readRDS(file = "//dapadfs/Workspace_cluster_9/Sustainable_Food_System/Results/sfs_index_knn_imputed.RDS")
indices[,1:(ncol(indices)-1)] <- round(indices[,1:(ncol(indices)-1)], 2)
n <- 4
colstops <- data.frame(
q = 0:n/n,
c = substring(viridis(n + 1), 0, 7)) %>%
list.parse2()
highchart(type = "map") %>%
hc_add_series_map(map = worldgeojson, df = indices, value = "SUFS", joinBy = "iso3") %>%
hc_colorAxis(stops = color_stops()) %>%
hc_tooltip(useHTML = TRUE, headerFormat = "",
pointFormat = "{point.name} has a SFS index of {point.SUFS}") %>%
hc_colorAxis(stops = colstops) %>%
hc_legend(valueDecimals = 0, valueSuffix = "%") %>%
hc_mapNavigation(enabled = TRUE) %>%
hc_add_theme(thm)
```
Column {.tabset data-width=400}
-----------------------------------------------------------------------
### PLS-PM 1 - Estimations
```{r}
complete_data <- readRDS(file = "//dapadfs/Workspace_cluster_9/Sustainable_Food_System/Input_data/data_joined.RDS")
# Method 1: k nearest neighbors (non-parametric alternative)
complete_data1 <- VIM::kNN(data = complete_data); complete_data1 <- complete_data1[,colnames(complete_data)]
# PLS-PM: Using repeated indicators
# Define path model matrix (inner model)
NUTR <- c(0, 0, 0, 0)
HINT <- c(0, 0, 0, 0)
FSCY <- c(0, 0, 0, 0)
SUFS <- c(1, 1, 1, 0)
sfs_path <- rbind(NUTR, HINT, FSCY, SUFS); rm(NUTR, HINT, FSCY, SUFS)
colnames(sfs_path) <- rownames(sfs_path)
sfs_blocks <- list(2:3, 4:5, 6:9, 2:9)
sfs_modes <- rep("A", 4)
sfs_pls <- plspm(complete_data1, sfs_path, sfs_blocks, modes = sfs_modes)
plot(sfs_pls)
```
### PLS-PM 1 - Indices correlations
```{r}
pairs(sfs_pls$scores, pch = 20)
```
### PLS-PM 1 - Crossloadings analysis
```{r}
xloads = melt(sfs_pls$crossloadings, id.vars = c("name", "block"))
gg <- ggplot(data = xloads, aes(x = name, y = value, fill = block))
gg <- gg + geom_hline(yintercept = 0, color = "gray75")
gg <- gg + geom_hline(yintercept = c(-0.5, 0.5), color = "gray70", linetype = 2)
gg <- gg + geom_bar(stat = 'identity', position = 'dodge')
gg <- gg + facet_wrap(block ~ variable)
gg <- gg + theme(axis.text.x = element_text(angle = 90), line = element_blank())
ggplotly(gg)
```
### Network analysis
```{r}
rownames(complete_data1) <- complete_data1$ISO3
# Calculate similarity measure
sfs_dis <- cluster::daisy(x = complete_data1[,-1], metric = c("gower"), stand = FALSE)
sfs_dis <- 1 - as.matrix(sfs_dis)
# Do cluster analysis
sfs_pca <- FactoMineR::PCA(X = complete_data1[,-1], scale.unit = T, graph = F)
sfs_hpc <- FactoMineR::HCPC(res = sfs_pca, nb.clust = -1, graph = F)
complete_data1$cluster <- sfs_hpc$data.clust$clust
# Visualize using networkD3
sfs_dis[lower.tri(sfs_dis, diag = TRUE)] <- NA
sfs_dis <- na.omit(data.frame(as.table(sfs_dis))); names(sfs_dis) <- c("from", "to", "similarity")
sfs_dis <- sfs_dis[sfs_dis$similarity >= .98,] # Filter by more than 98 degree of similarity
gD <- igraph::simplify(igraph::graph.data.frame(sfs_dis, directed = FALSE))
nodeList <- data.frame(id = c(0:(igraph::vcount(gD) - 1)), name = igraph::V(gD)$name) # because networkD3 library requires IDs to start at 0
getNodeID <- function(x){ which(x == igraph::V(gD)$name) - 1 } # to ensure that IDs start at 0
edgeList <- plyr::ddply(sfs_dis, .variables = c("from", "to", "similarity"),
function (x) data.frame(fromID = getNodeID(x$from),
toID = getNodeID(x$to)))
nodeList <- cbind(nodeList, nodeDegree = igraph::degree(gD, v = igraph::V(gD), mode = "all")); rm(gD, getNodeID)
nodeList$cluster <- as.numeric(as.character(complete_data1$cluster))[match(nodeList$name, complete_data1$ISO3)]
networkD3::forceNetwork(Links = edgeList,
Nodes = nodeList,
Source = "fromID",
Target = "toID",
Value = "similarity",
NodeID = "name",
Group = "cluster",
opacity = 1,
fontSize = 15)
```